Investigation of tracer diffusion in crowded cylindrical channel 
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Based on a coarse-grained model, we carry out molecular dynamics simulations to analyze the 
diffusion of a small tracer particle inside a cylindrical channel whose inner wall is covered with 
randomly grafted short polymeric chains. We observe an interesting transient subdiffusive behavior 
along the cylindrical axis at high attraction between the tracer and the chains, however, the long 
time diffusion is always normal. This process is found to be enhanced for the case that we immobilize 
the grafted chains, i.e. the sub-diffusive behavior sets in at an earlier time and spans over a longer 
time period before becoming diffusive. Even if the grafted chains are replaced with a frozen sea of 
repulsive, non-connected particles in the background, the transient subdiffusion is observed. The 
intermediate subdiffusive behavior only disappears when the grafted chains are replaced with a 
mobile background sea of mutually repulsive particles. Overall, the long time diffusion coefficient of 
the tracer along the cylinder axis decreases with the increase in system volume fraction, strength of 
attraction between the tracer and the background and also on freezing the background. We believe 
that the simple model presented here could be useful for a qualitative understanding of the process 
of macromolecular diffusion inside the nuclear pore complex. 
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I. INTRODUCTION 

Diffusion in crowded environment has been an active 
area of experimental and theoretical research p]4Tl| . Es- 
pecially in biology it is not uncommon to find aqueous 
crowded environments with agents that are important 
for biological functions. Many examples where crowding 
and/or sticky interaction between the diffusing species 
and the surroundings changes the rate and even the na- 
ture of the diffusion process have been found in biological 
or synthetic contexts: 

DNA-binding proteins search for specific binding sites 
on a DNA molecule by combining one dimensional diffu- 
sion along the DNA strand and three dimensional diffu- 
sion in the bulk [0, [l(| . The self diffusivity of proteins in 
crowded solutions has been shown to be slowed down to 
one fifth of the dilute limiting value . Also diffusion 
of proteins across the nuclear pore complex (NPC) (l2T - 
[l8j is an example of diffusion in a crowded environment 
where the crowding is due to the presence of proteins 
called nucleoporins that are rich in hydrophobic amino 
acids and form a brush [10, (2(| or a reversible hydrogel 
[H Ej. Proteins diffusing through NPC bind to these 
nucleoporins and experience a slowdown. A thorough 
explanation of this effect would lead to a better under- 
standing of the biological functioning of the NPC [8l.[2l|. 

Another issue is the emergence of non-Fickian, anoma- 
lous diffusion in presence of crowding Evidence of 
anomalous diffusion was found for colloidal tracer par- 
ticles in entangled actin filament networks Q and for 
soluble proteins in highly viscoelastic cytoplasm and nu- 
cleoplasm [22j. Very recently it has been shown that 
diffusion inside the NPC could also be anomalous [23| . 
On the other hand sticky interaction can also lead to 
subdiffusive behavior even in absence of crowding. It has 
been shown experimentally that functionalized colloidal 
particles with DNA sticky ends diffuse anomalously on 
a complementarily coated surface (24j . Diffusion in hy- 



drogels and gel-like media [3, 25 28] are also examples 
of diffusion in crowded environment and subdiffusive be- 
havior in an intermediate time scale has been observed in 
computer simulation [25| and signatures of subdiffusive 
behavior has been observed experimentally [3j. 

These observations suggest that similar physical princi- 
ples govern the diffusion processes in under very different 
conditions. We investigate these principles in a simple 
and generic model and shed light the role of crowding 
and sticky interaction with means of molecular dynamics 
simulations. 

We study the diffusion of tracer particles inside a cylin- 
drical channel grafted with polymeric chains from the in- 
side. The polymeric chains consists of mutually repulsive 
monomers connected with springs. In the field of polymer 
science, this is a typical model for a polymer in a good 
solvent [10] ■ Tracers are added to the system and the im- 
pact of their interaction with the polymer chain on their 
diffusional property is studied. The thermal effects of the 
surrounding water is replaced by a standard LAngevin 
thermostat. In this sense it is a solvent-implicit model 
where hydrodynamic effects are neglected. 

From a biological perspective this model, although be- 
ing very simple, is inspired by the brush model of the cen- 
tral plug of NPC [30| . Nucleoporins rich in hydrophobic 
units form a brush-like structure [2lj and it is believed 
that proteins being transported through the central plug 
bind to these nucleoporins [10, E[ • These bindings are 
believed to be hydrophobic in nature, and each binding 
is in the range of 1 — 2 ksT [32l. |33|. 

To mimic the aspect of crowding we chose a repul- 
sive interaction between all chain monomers and the 
tracer particles. To investigate the role of interaction 
we incorporate a Lennard- Jones potential between the 
tracer and the monomers of which we vary the interac- 
tion strength in a range from 1 — 2.5 fc^T corresponding 
to the interaction strength of hydrophobic contacts be- 
tween proteins and nucleoporins in the NPC. The result- 
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ing diffusive behaviour is characterized by investigating 
the mean-square-displacement for the tracer on different 
timescales. 

The paper is structured as follows. In section II, the 
model, the simulation methods and the classification of 
observed diffusion processes are discussed. In Section IV 
the results and interpretations for the different models 
are given. The paper ends with the conclusions in section 
IV. 



II. MODEL AND METHOD 

The model described in the following section is eval- 
uated using molecular dynamics. We use a generic unit 
system in which the Lennard- Jones parameter a is cho- 
sen as the unit of length, the temperature multiplied with 
Boltzmann's constant as the unit of energy and mass is 
measured in units of the mass of the monomers in the 
chain molecules. All quantities in the following are ex- 
pressed in this unit system. 

A cylinder with height 24<r and radius 9cr is created 
from particles with a diameter a and a mutual distance of 
1 er so that a closed cylindrical surface is formed on which 
each particle has four next neighbors. These particles, 
we call them wall particles, are fixed in space, i. e. their 
equations of motions are not integrated, independent of 
the force acting on them. The cylinder axis is chosen to 
be the z axis. Each of the polymeric chains is made of 
N = I 3 monomers, connected through a finite extensible 
nonlinear elastic (FENE) potential: 
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Here Kf is the force constant and r max is the maximum 
displacement of the bond. In the simulation the param- 
eters are chosen to be Kf = 7, r max = 2, N = 13. The 
cylindrical surface is grafted randomly with polymeric 
chains by fixing the end monomer of each of the chains 
to one of the wall particles. 

The repulsive interaction between the monomers them- 
selves and with the wall is modeled via the purely repul- 
sive variant of the Lennard- Jones potential also known 
as the Weeks-Chandler- Andersen (WCA) potential [H|]: 



WCA 



(r) = 




,if r < (2) 1 /V 
, otherwise, 



(2) 

with e = 1 and a = 1. 

The tracer particles are of the same size a and of the 
same mass as the monomers. They interact with the wall 
particles by means of the WCA potential, and with the 
grafted polymers either via the WCA potential (for sim- 
ple crowding), or by an attractive (= "usual") Lennard- 
Joncs potential: 
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, otherwise, 
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FIG. 1: Illustration of the simulation set up. Tracer parti- 
cles (white spheres) inside a cylindrical pore with rigid walls 
(black) grafted from inside with polymeric chains (blue). For 
graphical clarity the polymer chains are represented as (thin- 
ner) continuous curves (VMDs licorice drawing method) , and 
the size of the tracer particles is exaggerated. 



where the cut off radius is fixed at r cut = 2.5 and the 
attraction strength e is varied. 

The degree of crowding is characterized by the grafting 
density 7, i. e. the number of chains per unit are of the 
cylinder. A typical value of n = 75 chains corresponds 
to 7 = 0.055. 

A typical representation of the investigated system is 
shown in Fig. |T]by using VMD [35| . 

We performed five independent simulations for each 
data point with 2 million timesteps each in order to allow 
to estimate the statistical errors indicated in all figures. 
All simulations are carried out with the Langevin thermo- 
stat, thus solving the full phase space stochastic equation 
of motion known as the Langevin equation of motion. Pe- 
riodic boundary conditions along the cylindrical axis are 
used throughout the simulation and the velocity Vcrlct 
algorithm is used with a time step of r = 0.005. 

For a particle with mass m and the position co-ordinate 
r,; (where i stands for the zthe particle), interacting with 
all other particles (i ^ j) through a general poten- 
tial U (rij ) the full phase space (underdamped) Langevin 
equation of motion reads as 



dt 2 



D 



where 



Fr = -mr 



dcj(t) 
dt 



(4) 



(5) 



is the frictional drag acting on the particle with friction 
coefficient V, is the random force acting on the par- 
ticle. The first moment of the random force is chosen to 
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be zero: 

(Ff (i)> = 0. (6) 

The magnitude of the random forces obeys 

(Ff (ti)Ff (f a )> = erkBTmSySfa - f 2 ) (7) 

where kg is the Boltzmann constant and T is the tem- 
perature of the Langcvin thermostat. 

The delta function in time between the random forces 
in the above equation ensures the spectrum of the ran- 
dom force to be white, corresponding to a Markovian pro- 
cess in phase space with no memory. What is that: This 
also means dissipation is Ohmic. Throughout the simula- 
tion hydrodynamic interactions are neglected. With the 
above choice of random forces a fluctuation-dissipation 
theorems, that connects the diffusion coefficient D and 
the friction coefficient T holds: The diffusion coefficient 
of a single particle, up to a factor ksT/m is the inverse of 
the friction coefficient: Z?f rcc = ^tt = 1- In all simula- 
tions we apply a friction coefficient of 1 and the diffusion 
coefficients reported in the paper are given in units of 
Dfree, being also 1 in our unit system. 

The molecular dynamics simulation have been carried 
out with the help of the simulation package ESPResSo 
|36| in version 3.1 [13]. 

To describe the diffusion of tracer particles it is 
a common practice to use the mean square displace- 
ment (MSD), (Ar to (f) 2 ) as a measure, where Ar(f) = 
(f(t + to)— r(fo)) is the displacement of the particle at 
time t + t o from its initial position at io- The angular 
bracket denotes an average over all possible realizations 
of the process. Under the time translational invariance 
this quantity is independent of to and this index can be 
dropped, thus t can also be interpreted as the lag time. 
For ergodic systems, the average over realizations can be 
replaced by an ensemble average. In case of normal dif- 
fusion, MSD increases linearly with time and in d dimen- 
sion reads as (Ar(t) 2 ) = 2dDt, where D is the diffusion 
coefficient of the tracer particle. On the other hand, a 
process where the MSD is not linear in time but pro- 
portional to t a with a ^ 1 is referred to as anomalous 
or non-Fickian diffusion Here a is the diffusion ex- 
ponent and if a < 1, the process is called subdiffusive 
and if a > 1, the process is called supcrdiffusive. For 
anomalous diffusion there is no well defined diffusion co- 
efficient. While the definition given above is very helpful 
for isotropic environments, in anisotropic environments 
the diffusion constant needs to be replaced by a diffu- 
sion tensor. As we are only interested in the diffusive 
behaviour along the cylinder axis (=z-axis), we measure 
the MSD in z direction: 

MSD= ([r z (t + t Q )-r z (t Q )} 2 ), (8) 

corresponding to the zz component of the diffusion ten- 
sor. Throughout the paper D always represents the 
tracer diffusion coefficient along the cylinder axis and 
MSD is the mean square displacement of the tracer along 
the cylinder axis. 



III. RESULTS 

In this section, we present the results of the MD simu- 
lations. We compute the MSD of the tracer and analyze 
the resulting transient diffusive behaviour with respect 
to the question if anomalous diffusion appears and try to 
elucidate the reasons for this by comparison to reference 
systems. 

To explore the role of chain connectivity we compare 
all our results to a system when the chain monomers 
arc not connected via the FENE potential, but free to 
move thus forming a gas or a "sea" of independent parti- 
cles. The other reference is a system in which the grafted 
chains are moving much slower than the tracer parti- 
cle, we call this a frozen background. For brevity we 
characterize the four possible combinations by the key- 
words pairs (chains/particles) and (mobile/frozen). We 
will show that both properties, the chain connectivity 
and the background motion, have a distinct impact on 
the diffusion properties of the tracer particles. 

First we discuss the structure of the polymer coating: 
The radius of gyration of the chains is ~ 6 under all re- 
ported conditions and the mutual distance between graft- 
ing points is smaller also for the smallest grafting densi- 
ties. Thus we are in the regime where the grafted chains 
overlap considerably In Fig. [5] we report the local vol- 
ume fractions of chain monomers, i.e. the fraction of the 
volume locally occupied by monomers, assuming they are 
spheres of diamtcr a. The local volume fractions of up to 
10% and more also indicate that locally the crowding is 
significant: Under these conditions mean free paths ~ a 
are to be expected. 
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FIG. 2: The local monomer volume fraction, (j)(r) at differ- 
ent chain grafting densities (7) of the mobile grafted chains. 
Increasing the number of chains only affects the density but 
not the structure of the brush. 



A. Repulsive Tracer— Chain Interaction 

We now discuss the situation where the interaction be- 
tween tracer particles and grafted chains is purely repul- 
sive. We perform simulations with different numbers of 
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grafted chains (n) hence different grafting densities (7) 
and calculate the mean square displacement (MSD). We 
investigate four grafting densities seen in Fig. [3] In the 
double-logarithmic plot of the MSD, all lines nearly col- 
lapse (seen in the inset of Fig. [3]). We thus plot the MSD 
of the tracer along the cylinder axis divided by 2Dot, 
where Dq(= Df ree ) is the free diffusion coefficient along 
the cylinder axis. Plotting this reduced MSD allows to see 
the different characteristics of diffusion. On the timescale 
i a ballistic regime where the MSD is proportional to the 
square of t is observed, corresponding to a linear increase 
of the reduced MSD, and on longer timescales the diffu- 
sion becomes normal. As depicted in Fig. [3] the mean 
square displacement of the tracer particles shows no sign 
of anomalous diffusion for all investigated grafting den- 
sities. 
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FIG. 4: The relative diffusion coefficient D/Do for the purely 
repulsive tracers against the volume fraction of the system. 
The lines are guides to eyes. We compare the case where 
the monomers are connected by springs to form a chain to 
the case where they are isolated particles and the case where 
they are mobile vs. the case where they are frozen. 



infinite mass and is more likely to get bounced back to 
its initial position. At higher volume fractions the topol- 
ogy of the background becomes less important. As can 
be seen from the Fig. [4] the two data sets correspond- 
ing to frozen particles and chains coincide at higher vol- 
ume fractions as do the data sets for mobile particles and 
chains. 



Attractive Tracer-Chain Interaction 



FIG. 3: Comparison of MSD/2Dnf of the purely repulsive 
tracer along the cylinder axis of the cylinder with mobile 
grafted chains in the background at different chain grafting 
densities (7). The diffusion constant (seen as the t — s- 00 
limit of the curves) depends strongly on the density of grafted 
chains. 

The long time diffusion coefficient is reduced compared 
to a freely diffusing particle by collisions with the poly- 
mer chains. For low volume fractions of polymer chains 
we thus expect the diffusion coefficient to be reduced by 
a factor proportional to the (average) volume fraction, 
thus the probability to collide with a monomer per unit 
time. Neither the chain connectivity nor the question 
if the background moves on the same velocity scale as 
the tracer particles significantly influences this effect. In 
Fig. 2] we report the dependence of the diffusion coeffi- 
cient on the monomer volume fraction, d> = I" 1 ) Nn 
In all cases the tracer diffusion coefficient decreases sim- 
ilarly with increasing monomer volume fraction. If the 
monomers are connected and form chains, the impact 
on the diffusion coefficient is lower, as the volume from 
which the tracer is excluded is smaller. Also notice that 
the diffusion becomes faster with a mobile background 
than with a frozen background. This is because in the 
frozen background the tracer practically collides with an 



The dynamics of a tracer particle which interacts with 
the chains through an attractive Lennard- Jones potential 
as defined in (Eq.([3])) is more complex and shows inter- 
esting features that are missing in the case of a purely 
repulsive interaction. Very different modes of motion are 
important in case of attractive polymer-tracer interac- 
tion: An attractive tracer can get attached to one of 
the chains, follow the chains movement, move along the 
chain, hop from one chain to the other or hop back to 
the same chain. It can also move back to the bulk where 
it executes normal diffusion. The emerging process is a 
combination of all these modes of transport and should 
strongly depend on the strength of attraction between 
the tracer and the chain and also on the density of the 
grafted chains. We do not expect anomalous diffusion 
on long time scales, but see this rather as a transient 
phenomenon. This is what one would expect as long as 
the noise is white and the background is translationally 
invariant. In the simulations the strength of attraction 
between the tracer and the polymer (e) is varied from 1 to 
2.5 . The choice of this range of attraction is inspired by 
the fact that during the transport through NPC proteins 
bind to the nucleoporins through hy drop hobic contacts 
which are in the range of 1 — 2 fe^T. j33[. 

First we investigate the influence of the attractive in- 
teraction on the position where the tracer can be found. 
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In Fig. [5] the probability density to find the tracer at a 
given distance r from the cylinder axis for different at- 
traction strength at a fixed grafting density 0.055 is re- 
ported. Notice that the probability of finding the tracer 
along the radial direction (p(r)) changes profoundly on 
making the tracer-chain interaction attractive. In case 
of purely repulsive tracers, p(r) is high at the center of 
the cylinder (r = 0) and gradually decreases towards the 
chains and vanishes at the cylinder wall (r = 9). On the 
other hand in case of attractive tracers the density pro- 
file is qualitatively different. The tracer has a low density 
around the center of the cylinder (r = 0) and has a high 
density close to the chains. As expected the attractive 
interaction significantly pushes the tracers from the cen- 
ter of the cylinder towards the area of the grafted chains. 
For comparison we also insert the density profile 4>{r) of 
chain monomers. 
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FIG. 5: The probability of finding the tracer along the radial 
direction (p(r)) for different attraction strengths e at fixed 
grafting density 7 = 0.055. For comparison we repeat the 
density of monomers <j>(r) from Fig. [2] (solid line). For re- 
pulsive chain-tracer interactions the tracer is pushed towards 
the polymer-free center of the cylinder, while the attractive 
interaction pulls it into the polymer brush. 

Increasing the attraction strenght between the tracer 
particle and the chains in general slows down the dif- 
fusion. In table U we report the long time diffusion co- 
efficients obtained for different attraction strengths (e) 
at the grafting density 7 = 0.055. When the attraction 
strength is increased, a sub-diffusive regime for an inter- 
mediate time period is observed in the plot of the reduced 
MSD against t (Fig.© for e > 2k B T. 

We give the following interpretation: The energy land- 
scape seen by a tracer particle for a fixed configuration of 
the polymer chain has minima in which the particle gets 
trapped on these intermediate timescales. These minima 
are especially deep in vicinity of kinks in the polymer 
chain, where the attractive regions of several monomers 
overlap. This also occurs in vicinity of points, where two 
chains closely pass by each other. When a particle gets 
trapped, it moves together with the monomers respon- 
sible for the trapping. Since the chains are connected 
the corresponding displacement will only be small. The 



TABLE I: -jy of the tracer with grafted mobile chains in 
the background at the grafting density 7 = 0.055 for different 
attraction strengths (e). Increased attraction strength slows 
down diffusion. 
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0.38 ± 0.02 
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FIG. 6: Comparison of MSD/2Do£ along the cylinder axis 
for different cases for e = 2 at the grafting density 7 = 0.055. 
Subdiffusive transient behaviour indicated by a negative slope 
is seen in all cases, except for the case where the monomers 
are not connected and allowed to move. Chain structure and 
freezing the background facilitates subdiffusion. 



structure of the polymer chains itself is subject to Brow- 
nian motion, and this motion eventually helps the tracer 
to escape the traps. 

To justify this interpretation wc perform indepen- 
dent control experiments with a system where the chain 
monomers are not connected (a sea of particles) and in- 
vestigate the impact of strongly slowing down the motion 
of the background, hence look the motion of the tracer 
particle on a frozen background. We report here only the 
results for attraction strength e = 2, but qualitatively the 
results that were obtained for the other investigated at- 
traction strengths of 1, 1.5 and 2.5 kgT are similar. 

The observed reduced MSDs in all four cases are shown 
in Fig. [7] By giving up the chain structure of the environ- 
ment we expect to facilitate diffusion as the deep traps 
where the attractive regions of several monomers over- 
lap are less likely as entropy drives them away from each 
other. The long time diffusion coefficient is increased by 
a factor of 2.2 and the reduced MSD exhibits no subd- 
iffusive behaviour any more. Fixing the position of the 
disconnected monomers slows down diffusion again, as 
this suppresses the joint diffusion of complexes of tracers 
and monomers, but merely stops the tracers until they 
can escape from their traps. For e — 2 this stopping effect 



6 



is so pronounced that a transient sub-diffusive behaviour 
occurs (corresponding to the negative slope of the green 
curve in Fig. [7]), but for weaker attraction it disappears. 
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FIG. 7: Comparison of M SD /2D t along the cylinder axis for 
the purely repulsive case and different attraction strengths (e) 
with mobile grafted chains in the background. 

Freezing the structure of the chains does not change 
the distribution of traps, but their lifetimes. Trapped 
particles have to escape by their own Brownian motion 
and the random motion of the chains does no more fa- 
cilitate their escape due to random displacements of the 
traps. It is seen that a substantially more pronounced 
subdiffusivc regime occcurs. The long time diffusion co- 
efficient is decreased by 70 % by freezing the background. 
This observation qualitatively supports the proposition 
of Bickel and Bruinsma [H||. They showed with a very 
simple model for the transport across NPC that the fluc- 
tuating network of chains in the background acts as an 
extra noise and actually enhances the diffusion of a pro- 
tein immersed in it as compared to in a frozen network. 



observed on any timescale. If the tracer interaction, how- 
ever, is only weakly attractive (1-2.5 kgT) we observe a 
stronger slowdown than in the repulsive case, and also a 
pronounced subdiffusive regime. The subdiffusive regime 
is enlarged by freezing the background, corresponding to 
a much slower motion of the background particles. These 
findings arc summarized in table [XXJ 

Our findings suggest that the subdiffusive behaviour 
observed by Lowe et. al.[23| for the diffusion of tracers 
in the NPC is not caused by crowding from a polymer 
brush alone. Our results clearly demonstrate that simple 
attractive interactions between the crowding polymer- 
brush and the transported particles are sufficient to ob- 
serve transient subdiffusive behaviour. The presence of 
attractive interactions were recently revealed in the work 
of Lim et al. [39[ who could show that proteins that help 
transporting cargo through the NPC, do have an attrac- 
tive interaction with the nuclcoporins polymers in the 
NPC. It, however, remained unanswered if the nucleo- 
porins form a gel rather than a brush phase and what 
consequences this has for the cargo transport through 
the NPC. In the future we want to refine and extend our 
presented model into this direction. 



TABLE II: No Subdiffusion is observed in case of purely re- 
pulsive tracers. Results for attractive tracers are summarized 
below. 



Background Mobile 


Frozen 


Chains 


Transient 


Transient 




subdiffusion 


subdiffusion 




at intermedi- 


at intermedi- 




ate attraction 


ate attraction 


Particles 


No subdiffu- 


Transient 




sion 


subdiffu- 






sion at high 






attraction 



IV. CONCLUSIONS 

In this work, we have investigated the process of diffu- 
sion of tracer particles inside a crowded cylindrical chan- 
nel. Our coarse grained simulations show that if the in- 
teraction between the environment and the tracer parti- 
cle is purely repulsive, the diffusion is slowed down signif- 
icantly: Up to 80% if the background moves, up to 95% 
if the background is frozen. No subdiffusive behaviour is 
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